Pattern formation and selection in quasi-static fracture 
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Fracture in quasi-statically driven systems is studied by means of a discrete spring-block model. 
Developed from close comparison with desiccation experiments, it describes crack formation induced 
by friction on a substrate. The model produces cellular, hierarchical patterns of cracks, characterized 
by a mean fragment size linear in the layer thickness, in agreement with experiments. The selection 
of a stationary fragment size is explained by exploiting the correlations prior to cracking. A scaling 
behavior associated with the thickness and substrate coupling, derived and confirmed by simulations, 
suggests why patterns have similar morphology despite their disparity in scales. 

PACS numbers: 05.65.+b, 62.60.Mk, 46.50.+a 



Nature is full of fascinating patterns |Q . A ubiquitous 
yet relatively less explored class of patterns is that pro- 
duced by the fracture of solids M , as often seen in cracked 
structures, battered roads and dried out fields. Although 
there have been early observations and characterizations 
UHQ], only rather recently has it been studied system- 
atically The mathematical problem of a network 
of interacting cracks is a formidable one, in view of 
the difficulties facing a lone crack propagating in a ho- 
mogeneous medium [[l0| . One way to make progress is to 
turn to a statistical description. This view has been pur- 
sued extensively in investigations of the fracture of disor- 
dered media , borrowing concepts such as percolation 
and universality from studies of phase transitions. For 
crack patterns, an immediate observation is the geomet- 
rical similarities of patterns over a wide range of scales, 
from microns § to kilometers §. This suggests some 
universal mechanism is at work and microscopic details 
may be unimportant. Hence, analogous to phase tran- 
sitions, a mesoscopic, coarse-grained description may be 
sufficient and more useful than a microscopic one, pro- 
vided essential features are captured. 

Crack patterns often arise from slow physical (e.g., 
cmbrittlemcnt, contraction) or chemical (e.g., oxidation) 
variations, or their combination, in the material proper- 
ties of an overlayer on top of an inert substrate. The over- 
layer may fail in different ways [fl4j, such as decohesion, 
buckling, spalling and in-plane cracking. In this paper, 
we address the pattern-selection aspects of quasi-static, 
in-plane cracking by identifying the dominant mecha- 
nism and control parameters that determine the ensuing 
patterns. By understanding a simplified, coarse-grained 
model, we hope to achieve the same for the phenomena 
in general. 

A Bundle-Block Model - In a mesoscopic approach, 
the grains in the overlayer are represented by an array 
of blocks. Each pair of neighboring blocks is connected 
by a bundle of H bonds (coil springs) each of which 
has spring constant k = 1 and relaxed length I. Ini- 



tially, the blocks are randomly displaced by f = (x,y), 
where \r*\ <C I, about their mean positions on a triangular 
lattice. For slow cracking on a frictional substrate, the 
motion of the grains is overdamped, the system evolves 
quasi-statically with driving rate much slower than relax- 
ation rate. In this case, we need not solve the equations 
of motion but instead may update the configurations ac- 
cording to threshold criteria specified below. 

In typical crack formations, the layer often hardens 
and/or weakens in time; it tends to contract but is re- 
sisted by friction from the substrate. As a result, stress 
is built up and relaxed slowly. The simplest way to in- 
corporate these effects is first to pre-strain the array by 
s = (a — I) /a > 0, where a = 1 and I < a are the initial 
and final grain spacing, respectively. Then, two thresh- 
olds, F s for slipping and F c for cracking, are decreased 
systematically to induce slipping and cracking. This rule 
is more general than it seems: As the evolution of the 
model is determined entirely by the ratios between forces 
and thresholds, hardening, weakening and their combi- 
nation amount to this same rule. Since the thresholds de- 
crease by the same physical means, their ratio k = F c /F s 
remains fixed. In a drying experiment, the narrowing and 
breaking of liquid bridges between grains, due to evap- 
oration, may provide a physical picture of this driving. 
Even though there may not be an exact correspondence 
in experiments, and other drivings are conceivable [ p"3| |, 
our choice is physically relevant and sufficiently simple 
to allow for some analytic understanding. 

The general expression of the force is rather compli- 
cated, it results in slow updating. Since in realistic situ- 
ation the stress is primarily tensile and the strain s«l, 
to a good approximation we expand the force to first 
order in r to obtain a Hookean form for the force com- 
ponents on a block: 
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where the index i labels the six nearest neighbors and 
the associated bundles, and the constants Cj = cos[(2i — 
1)71-/6] and Si ~ sm[(2i — l)7r/6]. F y may be obtained 
from F x by interchanges x «-» y, and d <-> S^. 

Initially, the system is globally stable. In discrete sim- 
ulation steps t, it evolves according to these rules: 

1. The thresholds are lowered until cither is exceeded 
somewhere in the system; 



2. if F = ^F 2 + F 2 > F s , the block slips to a me- 
chanically equilibrium position where F = 0; 

3. if the bond tension > F c in a bundle of Hi bonds, 
Hi — > Hi — 1 provided the bundle has not been 
damaged in the same step. 

The system relaxes until global stability is restored. 
This constitutes one simulation step or one event. In 
rule 3, we break no more than one bond/bundle/step 
in order to avoid instantaneous downward propagation 
of cracks before the neighboring stress field has ever re- 
laxed. Otherwise, spurious long and vertical microcracks 
would be generated which do not conform to real systems. 
The progressive damage of the bundle mimics successive 
breakage of liquid bridges, hence H plays the role of the 
local thickness . Our algorithm is designed to take ad- 
vantage of the dominance of the horizontal propagation 
of cracks, due to the presence of the top, free surface. 

Results and Discussions - Throughout this work, 
square systems are used with linear size 20 < L < 300, 
s < 0.1, and free boundary conditions (FBC) are im- 
posed p6| . In the k — H subspacc, the time evolution 
can be divided into three qualitatively distinct phases: 

I. The system contracts by slippings, the slip size (to- 
tal number of slippings in an event) grows in time; 

II. There are slippings and bond breakings, the system 
is progressively damaged, then fragmented, while 
contraction continues; 

III. Bond-breaking saturates, fragmentation stops, slip- 
ping dominates as block spacing converges to the 
equilibrium value I. 

A typical evolution is shown in Fig. [l]. The main fea- 
ture of phase I is the growth of correlations. Since ev- 
erything depends on F s , we need to know how it varies 
first. Its decay rate is governed by the density of F near 
F s among the L 2 blocks, where < F < F s . To a first 
(mean-field) approximation, this density is F s /L . Hence 
F s (t) -F s (t + 1) cx F s (t)/L 2 , leading to F s {t) cx e~^ T 
with a decay constant r ~ L z and a mean-field prediction 
zmf — 2. Numerically, we find z = 1.80 ± 0.03, reflect- 
ing the finite, albeit weak, spatial correlation among the 
seeds that initiate distinct events. 



Now we derive the growth laws for the strain field and 
the slip size. Due to the diminishing threshold, slippings 
start from the free edge and invade into the bulk, giving 
rise to a strain-relieved peripheral region of width de- 
noted by £. The strain increases from the boundary to- 
ward the bulk, reaching the bulk value s (see Fig. H(a)). 
To determine its profile, consider for the moment a one- 
dimensional version: the force on the i th block is given 
by Fi — (a — I + — Xi)H — (a — I + Xi — Xi-±)H = 
(sj — Sj_i)if. Hence F = Hds/du, where u is the dis- 
tance from the edge. Metastability of the system requires 
Hds/du < F s but of the same order. Thus we obtain for 
general dimensions the growth of the invasion length 
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with constant £0 = 0(1). For a finite system of size L, 
the invasion saturates (i.e., £ = L) when F s < F s sat = 
£ sH/L, or equivalently t > t sat = r ln(i//F| at ). As a 
result, the global maximum strain s max in the system 
drops continuously thereafter (see Fig. 0(b)): 
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As grows, so does the mean slip size, S(t), as re- 
quired by energy balance. The total energy relieved over 
the strain-relieved region is of the order U = ^Ls 2 H/2. It 
is dissipated by friction during slippings. Since the slip 
distance is aF s /H || where a = 1/3(1 + s) Q, each 
block slip dissipates an energy aF 2 /2H. The total dissi- 
pation is then E = J* dt' S{t')aF s (t') 2 /2H . Equating E 
to U and differentiating with respect to t, we readily find 



S(t) = So 



H 

m 



t < t s 



(4) 



where the prefactor So — £ l os 3 L/2ar. Since the invasion 
is independent of L, we have £0 ~ L° and hence Sq ~ 
L 1_z . For t > t S at, S(t) saturates to a constant S sat = 
So(L/st;o) 3 ~ L n which diverges as L — > 00: the system 
becomes critical with a nontrivial exponent 77 > 0. Put 
together, we obtain a "scaling relation" 
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relating the dynamical and stationary behavior. From 
zmf = 2, we get ?7mf = 2. Simulations give r/ = 
2.20 ± 0.03 instead, in agreement with the independent 
estimate of z. From Eq. (||), it is evident that dynamic 
scaling S(t;H,L) = U><f>(H/LF s (t)) is obeyed, with the 
scaling function &(x <C 1) ~ x 3 and $(x ^> 1) =const, 
as exhibited in Fig. [| 

Several remarks are in order: (i) 77 > 2 implies that 
each site slips an average L v ~ 2 times per event. Such re- 
peated topplings are manifestations of system- wide corre- 
lations in the force variable, (ii) The connection between 
dynamics (z) and criticality (77) is attributed to energy 
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balance. Thus we believe similar relations also exist for 
other forms of driving, (iii) For finite k, fracture eventu- 
ally sets in before criticality is fully developed. However, 
even then, each fragment settles into a stationary, critical 
state characterized by the same rj and its fragment size. 

Moving on to phase II, an important question in frac- 
ture is how its onset is affected by material properties and 
external conditions. Apparently, fracture occurs when 
kF s decreases below the bond tension s in the bulk, i.e., 

F s -^F;^- (6) 

K 

from above. Note that this result is independent of 
H, but the corresponding simulation step is: t* ~ 
rln(reff/s), because F s w He~ t / T . Hence, the substrate 
and the thickness have the same effect on the onset. The 
extra waiting time for thicker layers or more frictional 
substrates allows more slippings and hence stronger cor- 
relations at t* - larger fragments are anticipated. 

The difference in correlations is also reflected in the 
morphology of cracks. For very small reiJ, the substrate 
is very frictional or the layer very thin. Pinning is so 
strong that most avalanches are of size one. In a virtu- 
ally uncorrelated stress environment, cracks soon appear 
but do not propagate. Fragmentation is resulted from 
percolation of microcracks. For kH <C 1, most bonds are 
ultimately broken to turn the system into powder. 

Going up in kH, the larger friction or thickness allows 
stronger correlations up to fracture, measured by 
and S(t*). Relaxation of stress field around crack tips be- 
comes more effective, leading to more pronounced stress 
concentration and straight crack propagation. The re- 
sulting fragmentation process is hierarchical: after a cel- 
lular network is formed with the primary straight cracks, 
secondary cracks break the fragments into smaller ones, 
and so on. Cracks of later generations are more "diffu- 
sive" in character, they tend to wiggle and branch under 
smaller stresses Jl8|, as there are more and more free 
boundaries. These features are exemplified in Fig. || At 
long times, fragmentation is complete, resulting in frag- 
ments with areas that fit a log-normal distribution (cf. 
||), and a mean area A ~ H 2 (see below). Experimen- 
tally, morphology may be characterized by the distribu- 
tion P(9) of crack angle 9 at a joint. A shift of the 
peak position from 9 « 90° to 60° as thickness decreases 
has been observed in amorphous media Q . Although the 
usefulness of P(9) is limited by lattice anisotropics in our 
model, a direct comparison with real patterns (see Fig. |J) 
produced by slow drying of cornstarch-water mixtures (to 
be reported elsewhere) reveals striking similarities. Fur- 
thermore, the linear dependence of A on thickness agrees 
with previous experiments [BR . 

As kH is increased further, the question of fracture or 
not arises. Fracture is possible only if nF s (t) < s max (f) at 
some t. From Fig. |^(b), this translates to the condition 
nF s (t sat ) = s (cquivalently t sat = t*). A critical system 
size L c (k, H) = £,qkH emerges, below which a system 



never cracks at fixed re and H. Equivalently, a finite sys- 
tem of size L does not crack if kH > L/£q. This again 
agrees with the basic experimental fact that a sufficiently 
small and thick layer does not crack. Now, given that a 
system of size L > L c cracks, the next, even more inter- 
esting question of pattern selection is when fragmentation 
stops. We do not have a rigorous answer, as the dynam- 
ics of interacting cracks is complicated, but we argue that 
it stops when the system has broken up into fragments 
of typical size 0{L C ). This gives A(k,H) - (reiJ) 2 . This 
argument implies the memory of the correlations at the 
onset persists in the final fragmented state, reminiscent 
of the generation of quenched disorder in systems with 
metastability and frustration [jl9| - both features crucial 
to surface fracture. Excellent scaling behavior in simula- 
tions (Fig. ||) strongly supports this scenario. 

In summary, we study a bundle-block model for quasi- 
static, in-plane fracture. Guided by experiments, we in- 
corporate two factors that dictate the pattern-forming 
process - the substrate and thickness - by means of two 
control parameters. The ensuing pattern, however, is 
characterized mostly by one scaling variable. Self orga- 
nization of the strain field prior to fracture is derived and 
shown to dictate the properties of the pattern, particu- 
larly to select the fragment size. Our results answer, at 
least qualitatively, why thin layers crack finely and thick 
ones do not crack at all, and why similarities are accom- 
panied by a variability of scales. Supported by experi- 
mental observations, our model appears to capture the 
salient features of real systems, so that despite the use of 
a specific form of driving, the ideas may apply to other 
situations generally dominated by stick-slip motion. 
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FIG. 1. Snapshots for L = 100, s = 0.1, k = 0.5, H = 9. 
Note the progression from straight to diffusive cracks. 
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FIG. 2. (a) Strain field profiles show progressive strain 
relief as t increases for L = 100, s = 0.1, H = 2, k = oo. It 
leads to temporal variation of the maximum strain in (b). 
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FIG. 3. Dynamic finite-size scaling of mean slip size to- 
ward criticality. s — 0.1, k = oo, rj = 2.2 ± 0.03. 
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FIG. 4. Stationary patterns at different thickness. Up- 
per row: by simulation, L = 100, s = 0.1, k = 0.25. 
Increasing k at fixed H has the same effect. Lower row: 
from slow-drying experiment of cornstarch-water mixture 
with thickness matching that of simulation. 
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